bulk RNAseq basics
Jose Alejandro Romero Herrera
2021-11-18
1 Introduction
This is a basic bulk RNAseq analysis using the package DESeq2 and gprofiler2 for Differential Expression Analysis and Functional Annotation. The code is a modified version of the DESeq2 and gprofiler2 vignettes. The data used for this example below is the pasilla dataset: an RNAseq experiment that studied the effect of RNAi knockdown of Pasilla, the Drosophila melanogaster ortholog of mammalian NOVA1 and NOVA2, on the transcriptome.
2 Load libraries
library(tidyverse, quietly = T)
library(DESeq2, quietly = T)
library(RColorBrewer, quietly = T)
library(pheatmap, quietly = T)
library(gprofiler2, quietly = T)
library(Glimma, quietly = T)3 Load data
In order to proceed with the analysis we need two dataframes: - Count matrix dataframe with gene names as rows (and rownames) and samples as columns (and column names). - Sample metadata dataframe with samples as rows (and rownames) and metadata variables as columns (and column names).
Additionally, you can add information about your features (genes), with genes as rows and gene metadata on columns.
Ideally, the sample metadata will contain information about the conditions, treatments, replicates, etc.
The pasilla dataset contains a count matrix and the metadata, which is a good example to use as a template. The count matrix looks like this:
cts <- data.frame(read.csv(file = "../../data/example_data/pasilla_cts.tsv", sep="\t", row.names="flybase_id"))
head(cts)While the sample metadata looks like this:
coldata <- read.csv(file = "../../data/example_data/pasilla_metadata.tsv", row.names=1, sep = "\t")
head(coldata)Making the metadata factors will establish the order of the metadata variables. If you never tell the DESeq2 functions which level you want to compare against (e.g. which level represents the control group), the comparisons will be based on the alphabetical order of the levels. This is important since the first level of the metadata will be the reference sample/condition for the DE analysis.
coldata$condition <- factor(coldata$condition)
paste("The order of the conditions is:", paste(unique(coldata$condition), collapse=", "))## [1] "The order of the conditions is: untreated, treated"coldata$type <- factor(coldata$type)
paste("The order of the types is:", paste(unique(coldata$type), collapse=", "))## [1] "The order of the types is: single-read, paired-end"3.1 Metadata factor levels
If you want to change the order of the factor levels you can use one of these two. NOTE: this has to be done before running the DESeq() function.
#coldata$condition <- factor(coldata$condition, levels = c("untreated","treated"))
coldata$condition <- relevel(coldata$condition, ref = "untreated") #Specifies "untreated" as the reference level3.2 Create DESeq object
We create the DESeq object from the count matrix and the metadata. We specify that the analysis will be based on the design column.
dds <- DESeqDataSetFromMatrix(countData = cts,
colData = coldata, rowData = rownames(cts),
design = ~ condition)
dds## class: DESeqDataSet
## dim: 14599 7
## metadata(1): version
## assays(1): counts
## rownames(14599): FBgn0000003 FBgn0000008 ... FBgn0261574 FBgn0261575
## rowData names(1): X
## colnames(7): untreated1 untreated2 ... treated2 treated3
## colData names(2): condition type3.3 Pre-filtering
We can reduce the computational resources of the analysis by removing genes that are very hardly expressed. Additionally, you can collapse technical replicates using the function collapseReplicates().
keep <- rowSums(counts(dds)) >= 10
dds <- dds[keep,]3.4 Normalization
Finally the DESeq() function will normalize your read counts and estimate size factors per genes. This function is essential for the rest of the analysis.
dds <- DESeq(dds)4 Exploratory analysis
After running the DESeq() function, we can start our RNAseq analysis by exploring the relationships between our samples. In principle, all replicates for a specific condition should be similar to each other after normalization. Should this not be the case, one might need to remove outlier replicates. This is one of the reasons why RNAseq experiments need at least 3 replicates per condition.
4.1 Data transformations
For visualization and clustering (exploratory analyses) – it might be useful to work with transformed versions of the count data. There are two main methods used for this purpose: variance stabilizing transformations (VST), and the regularized logarithm or rlog. Both transformations produce transformed data on the log2 scale which has been normalized with respect to library size or other normalization factors.
vsd <- vst(dds) #Variance Stabilizing Transformation, vst is faster with larger number of samples
rld <- rlog(dds) #Regularized 4.2 Samples comparisons
We can do some sanity checks of the samples. We can see how they correlate to each other in a heatmap or see their Principal Components with a PCA plot.
4.2.1 Heatmap of the sample-to-sample distances
This plot shows how far away are each sample from each other. The darker the blue, the closer they are.
sampleDists <- dist(t(assay(vsd)))
sampleDistMatrix <- as.matrix(sampleDists)
rownames(sampleDistMatrix) <- paste(vsd$condition, vsd$type, sep="-")
colnames(sampleDistMatrix) <- NULL
colors <- colorRampPalette( rev(brewer.pal(9, "Blues")) )(255)
pheatmap(sampleDistMatrix,
clustering_distance_rows=sampleDists,
clustering_distance_cols=sampleDists,
col=colors)
4.2.2 Principal component plot of the samples
PCA plot using the first two components
pcaData <- plotPCA(vsd, intgroup=c("condition", "type"), returnData=TRUE)
percentVar <- round(100 * attr(pcaData, "percentVar"))
ggplot(pcaData, aes(PC1, PC2, color=condition, shape=type)) +
geom_point(size=3) +
xlab(paste0("PC1: ",percentVar[1],"% variance")) +
ylab(paste0("PC2: ",percentVar[2],"% variance")) +
coord_fixed() + theme_bw()
### Glimma plots Interactive visualizations of the DESeq results using the Glimma package, which provides excellent options for MA, Volcano and dimensionality reduction (like a PCA) plots. Use the groups argument to provide the condition or factor of your experiment. Unfortunately, the interactive plots created here are not fully compatible when knitting this document. Feel free to explore them using this link!
#glimmaMDS(dds, groups = dds$condition)Save your interactive plots using the htmlwidgets packages, specifically, the saveWidget()
#htmlwidgets::saveWidget(glimmaMDS(dds, groups = dds$condition), "mds-plot.html")4.3 Plot counts
It can be useful to examine the counts of reads for a single gene across the conditions. A simple function for making this plot is plotCounts, which normalizes counts by the estimated size factors (or normalization factors if these were used) and adds a pseudocount of 1/2 to allow for log scale plotting. The counts are grouped by the variables in intgroup argument, where more than one variable can be specified. You can select the gene to plot by its name or by numeric index.
d <- plotCounts(dds, gene= 1, intgroup="condition",
returnData=TRUE)
ggplot(d, aes(x=condition, y=count, color = condition)) +
geom_point(position=position_jitter(w=0.1,h=0)) +
scale_y_log10(breaks=c(25,100,400)) + theme_bw()
5 Differential expression analysis
In order to test for differential expression, we operate on raw counts and use discrete distributions After checking that our samples are well correlated within conditions, we can proceed to test for Differantially Expressed Genes using the results() function.
res <- results(dds)
head(data.frame(res))Information about the results columns can be retrieved using the following snippet:
mcols(res)$description## [1] "mean of normalized counts for all samples"
## [2] "log2 fold change (MLE): condition treated vs untreated"
## [3] "standard error: condition treated vs untreated"
## [4] "Wald statistic: condition treated vs untreated"
## [5] "Wald test p-value: condition treated vs untreated"
## [6] "BH adjusted p-values"In this example, we only have two conditions: “untreated” and “treated”, so there is no need to specify which comparison we would like to make. If your experiment contains more than two conditions, there are different ways of specifying which comparison you would like to make.
First, the “contrast” argument will take a vector that specifies the metadata and the order of the comparison. In this case, you can swap the order of the condition levels.
res <- results(dds, contrast=c("condition","untreated","treated"))
head(data.frame(res))Second, “name” argument will use the results of the resultsNames() function. This function creates possible combinations of conditions, but always using the reference condition for the comparison.
resultsNames(dds)## [1] "Intercept" "condition_treated_vs_untreated"res <- results(dds, name="condition_treated_vs_untreated")
head(data.frame(res))The notebook contrast_design.Rmd contains information about how to understand and define contrasts using model matrices, which might be a more intuitive way of comparing samples or conditions.
5.1 Extracting results
The output of the res() function is the one you want to use/save to identify differentially expressed genes using log2 Fold Change or adjusted p-value thresholds.
There is a chance that the adjusted p-value is NA. If you will work on this results, it might be helpful to change all NA to an adjusted p-value of 1.
res$padj[is.na(res$padj)] <- 1We can have a quick look of the adjusted p-value results by using the following custom function:
results_summary <- function(x, alpha = 0.05, LFC = 1) {
ngenes <- nrow(x)
signif <- sum(x$padj < alpha, na.rm = T)
up <- sum(x$padj < alpha & x$log2FoldChange > LFC, na.rm = T)
down <- sum(x$padj < alpha & x$log2FoldChange < -LFC, na.rm = T)
results <- c(paste0("Number of genes: ", ngenes),
paste0("Number of genes with adjusted p-value < ",alpha,": ", signif, " (", round((signif/ngenes)*100,digits = 2),"%)"),
" Of those:",
paste0(" with LFC < ", -LFC, ": ", down, " (", round((down/ngenes)*100,digits = 2),"%)"),
paste0(" with LFC > ", LFC, ": ", up, " (", round((up/ngenes)*100,digits = 2),"%)"))
writeLines(paste(results, collapse = "\n"))
return(paste(results, collapse = "\n"))
}res_summary <- results_summary(res, alpha = 0.05, LFC = 1)## Number of genes: 9921
## Number of genes with adjusted p-value < 0.05: 841 (8.48%)
## Of those:
## with LFC < -1: 107 (1.08%)
## with LFC > 1: 115 (1.16%)5.1.1 Save significant results
Extracting results based on adjusted p-value and LFC. this way we will get a filteres table with genes with LFC > 1 or LFC < -1 with an adjusted p-value of 0.05
LFC <- 1
adj_pvalue <- 0.05
sig_res <- res[res$padj < adj_pvalue & abs(res$log2FoldChange) > LFC,]
#write.table(x = sig_res, "./significant_results.tsv", quote = F, col.names = T, row.names = T, sep = "\t")6 DEA visualizations
6.1 Count matrix heatmaps
We can visualize different heatmaps depending on what we are interested on. We can see the top most expressed genes, the top most variable genes and the DE genes.
6.1.1 Top most variable genes
ntop <- 20
select <- order(rowSds(counts(dds,normalized=TRUE)),
decreasing=TRUE)[1:ntop]
df <- as.data.frame(colData(dds)[,c("condition","type")])
pheatmap(assay(vsd)[select,], cluster_rows=T, show_rownames=FALSE, scale = "row",
cluster_cols=FALSE, annotation_col=df)
6.1.2 Top expressed genes
ntop <- 20
select <- order(rowMeans(counts(dds,normalized=TRUE)),
decreasing=TRUE)[1:ntop]
df <- as.data.frame(colData(dds)[,c("condition","type")])
pheatmap(assay(vsd)[select,], cluster_rows=FALSE, show_rownames=FALSE,
cluster_cols=FALSE, annotation_col=df)
6.1.3 Differentially expressed genes
select <- rownames(sig_res)
df <- as.data.frame(colData(dds)[,c("condition","type")])
pheatmap(assay(vsd)[rownames(vsd) %in% select,], cluster_rows=TRUE, show_rownames=FALSE, scale = "row",
cluster_cols=FALSE, annotation_col=df)
6.2 Fold change plots
For fold change plots, it is useful to shrunk log2 fold changes, which removes the noise associated with log2 fold changes from low count genes without requiring arbitrary filtering thresholds.
resLFC <- lfcShrink(dds, coef="condition_treated_vs_untreated", type="apeglm")
resLFC$padj[is.na(resLFC$padj)] <- 16.2.1 MA plot
alpha <- 0.05
ggplot(data.frame(resLFC)) + theme_bw() +
labs(title = "MA plot", subtitle = res_summary,
color = paste0("Adjusted p-value < ",alpha)) +
geom_point(aes(x = log10(baseMean), y = log2FoldChange,
color = factor((padj < alpha), levels = c("TRUE","FALSE"))), size = 1)
6.2.2 Volcano plot
alpha = 0.05
LFC = 1
ggplot(data.frame(resLFC)) + theme_bw() +
labs(color ="Significant genes", title = "Volcano plot", subtitle = res_summary) +
geom_hline(yintercept = -log10(alpha), linetype = "dashed") +
geom_vline(xintercept = c(-LFC,LFC), linetype = "dashed") +
geom_point(aes(x = log2FoldChange, y = -log10(padj),
color = factor((padj < alpha & abs(log2FoldChange) > LFC), levels = c("TRUE","FALSE"))), size = 1)
6.2.3 Glimma plots
We can visualize our results interactively using the glimmaMA() and glimmaVolcano() functions.
#glimmaMA(dds, groups = dds$condition)
#glimmaVolcano(dds, groups = dds$condition)Save your interactive plots using the htmlwidgets::saveWidget() function.
#htmlwidgets::saveWidget(glimmaMA(dds, groups = dds$condition), "ma-plot.html")
#htmlwidgets::saveWidget(glimmaVolcano(dds, groups = dds$condition), "volcano-plot.html")7 Functional analysis with gprofiler2
gost() function allows us to do functional profiling of gene lists, such as our differentially expressed genes. The function performs statistical enrichment analysis to find over-representation of terms from Gene Ontology, biological pathways like KEGG and Reactome, human disease annotations, etc. This is done by using hypergeometric tests that are corrected for multiple testing.
7.1 Single query
A standard input of the gost() function is a (named) list of gene identifiers. The list can consist of mixed types of identifiers (proteins, transcripts, microarray IDs, etc), SNP IDs, chromosomal intervals or functional term IDs.
The result is a named list where result is a data.frame with the enrichment analysis results and meta containing a named list with all the metadata for the query.
gostres <- gost(query = rownames(sig_res),
organism = "dmelanogaster", ordered_query = FALSE,
multi_query = FALSE, significant = FALSE, exclude_iea = FALSE,
measure_underrepresentation = FALSE, evcodes = FALSE,
user_threshold = 0.05, correction_method = "g_SCS",
domain_scope = "annotated", custom_bg = NULL,
numeric_ns = "", sources = NULL, as_short_link = FALSE)
head(gostres$result)The result data.frame contains the following columns:
names(gostres$meta)## [1] "query_metadata" "result_metadata" "genes_metadata" "timestamp"
## [5] "version"7.2 Multiple queries
The function gost() also allows to perform enrichment on multiple input gene lists. Multiple queries are automatically detected if the input query is a list of vectors with gene identifiers and the results are combined into identical data.frame as in case of single query.
multi_gostres1 <- gost(query = list("chromX" = c("X:1000:1000000", "rs17396340",
"GO:0005005", "ENSG00000156103", "NLRP1"),
"chromY" = c("Y:1:10000000", "rs17396340",
"GO:0005005", "ENSG00000156103", "NLRP1")),
multi_query = FALSE)
head(multi_gostres1$result, 3)The column “query” in the result dataframe will now contain the corresponding name for the query. If no name is specified, then the query name is defined as the order of query with the prefix “query_.” Another option for multiple gene lists is setting the parameter multiquery = TRUE. Then the results from all of the input queries are grouped according to term IDs for better comparison.
multi_gostres2 <- gost(query = list("chromX" = c("X:1000:1000000", "rs17396340",
"GO:0005005", "ENSG00000156103", "NLRP1"),
"chromY" = c("Y:1:10000000", "rs17396340",
"GO:0005005", "ENSG00000156103", "NLRP1")),
multi_query = TRUE)
head(multi_gostres2$result, 3)7.3 Visualization
The enrichment results are visualized with a Manhattan-like-plot using the function gostplot() and the previously found gost results gostres:
gostplot(gostres, capped = TRUE, interactive = TRUE)The function publish_gostplot() takes the static plot object as an input and enables to highlight a selection of interesting terms from the results with numbers and table of results. These can be set with parameter highlight_terms listing the term IDs in a vector or as a data.frame with column “term_id” such as a subset of the result dataframe.
First we create the static plot
p <- gostplot(gostres, capped = FALSE, interactive = FALSE)
p
Then we make it in high quality. We can add highlighted terms if we want with the highlight_terms argument.
pp <- publish_gostplot(p, highlight_terms = c("GO:0005372"),
width = NA, height = NA, filename = NULL )
The gost results can also be visualized with a table. The publish_gosttable() function will create a nice-looking table with the result statistics for the highlight_terms from the result data.frame. The highlight_terms can be a vector of term IDs or a subset of the results.
publish_gosttable(gostres, highlight_terms = gostres$result[c(1:10),],
use_colors = TRUE,
show_columns = c("source", "term_name", "term_size", "intersection_size"),
filename = NULL) 
8 Session info
Finally, we create a session_info() table that will allow anyone to check what versions of R and packages are we using for reproducibility purposes.
devtools::session_info()## ─ Session info 💪🏾 👗 🦵🏽 ─────────────────────────────────────────────────
## hash: flexed biceps: medium-dark skin tone, dress, leg: medium skin tone
##
## setting value
## version R version 4.1.0 (2021-05-18)
## os macOS Big Sur 10.16
## system x86_64, darwin17.0
## ui X11
## language (EN)
## collate en_US.UTF-8
## ctype en_US.UTF-8
## tz Europe/Copenhagen
## date 2021-11-18
## pandoc 2.14.0.1 @ /usr/local/bin/ (via rmarkdown)
##
## ─ Packages ───────────────────────────────────────────────────────────────────
## package * version date (UTC) lib source
## annotate 1.72.0 2021-10-26 [1] Bioconductor
## AnnotationDbi 1.56.2 2021-11-09 [1] Bioconductor
## apeglm 1.16.0 2021-10-26 [1] Bioconductor
## assertthat 0.2.1 2019-03-21 [1] CRAN (R 4.1.0)
## backports 1.3.0 2021-10-27 [1] CRAN (R 4.1.0)
## bbmle 1.0.24 2021-08-05 [1] CRAN (R 4.1.0)
## bdsmatrix 1.3-4 2020-01-13 [1] CRAN (R 4.1.0)
## Biobase * 2.54.0 2021-10-26 [1] Bioconductor
## BiocGenerics * 0.40.0 2021-10-26 [1] Bioconductor
## BiocParallel 1.28.0 2021-10-26 [1] Bioconductor
## Biostrings 2.62.0 2021-10-26 [1] Bioconductor
## bit 4.0.4 2020-08-04 [1] CRAN (R 4.1.0)
## bit64 4.0.5 2020-08-30 [1] CRAN (R 4.1.0)
## bitops 1.0-7 2021-04-24 [1] CRAN (R 4.1.0)
## blob 1.2.2 2021-07-23 [1] CRAN (R 4.1.0)
## broom 0.7.10 2021-10-31 [1] CRAN (R 4.1.0)
## bslib 0.3.1 2021-10-06 [1] CRAN (R 4.1.0)
## cachem 1.0.6 2021-08-19 [1] CRAN (R 4.1.0)
## callr 3.7.0 2021-04-20 [1] CRAN (R 4.1.0)
## cellranger 1.1.0 2016-07-27 [1] CRAN (R 4.1.0)
## cli 3.1.0 2021-10-27 [1] CRAN (R 4.1.0)
## coda 0.19-4 2020-09-30 [1] CRAN (R 4.1.0)
## codetools 0.2-18 2020-11-04 [1] CRAN (R 4.1.0)
## colorspace 2.0-2 2021-06-24 [1] CRAN (R 4.1.0)
## crayon 1.4.2 2021-10-29 [1] CRAN (R 4.1.0)
## crosstalk 1.2.0 2021-11-04 [1] CRAN (R 4.1.0)
## data.table 1.14.2 2021-09-27 [1] CRAN (R 4.1.0)
## DBI 1.1.1 2021-01-15 [1] CRAN (R 4.1.0)
## dbplyr 2.1.1 2021-04-06 [1] CRAN (R 4.1.0)
## DelayedArray 0.20.0 2021-10-26 [1] Bioconductor
## desc 1.4.0 2021-09-28 [1] CRAN (R 4.1.0)
## DESeq2 * 1.34.0 2021-10-26 [1] Bioconductor
## devtools 2.4.2 2021-06-07 [1] CRAN (R 4.1.0)
## digest 0.6.28 2021-09-23 [1] CRAN (R 4.1.0)
## dplyr * 1.0.7 2021-06-18 [1] CRAN (R 4.1.0)
## edgeR 3.36.0 2021-10-26 [1] Bioconductor
## ellipsis 0.3.2 2021-04-29 [1] CRAN (R 4.1.0)
## emdbook 1.3.12 2020-02-19 [1] CRAN (R 4.1.0)
## evaluate 0.14 2019-05-28 [1] CRAN (R 4.1.0)
## fansi 0.5.0 2021-05-25 [1] CRAN (R 4.1.0)
## farver 2.1.0 2021-02-28 [1] CRAN (R 4.1.0)
## fastmap 1.1.0 2021-01-25 [1] CRAN (R 4.1.0)
## forcats * 0.5.1 2021-01-27 [1] CRAN (R 4.1.0)
## fs 1.5.0 2020-07-31 [1] CRAN (R 4.1.0)
## genefilter 1.76.0 2021-10-26 [1] Bioconductor
## geneplotter 1.72.0 2021-10-26 [1] Bioconductor
## generics 0.1.1 2021-10-25 [1] CRAN (R 4.1.0)
## GenomeInfoDb * 1.30.0 2021-10-26 [1] Bioconductor
## GenomeInfoDbData 1.2.7 2021-11-16 [1] Bioconductor
## GenomicRanges * 1.46.0 2021-10-26 [1] Bioconductor
## ggplot2 * 3.3.5 2021-06-25 [1] CRAN (R 4.1.0)
## Glimma * 2.4.0 2021-10-26 [1] Bioconductor
## glue 1.5.0 2021-11-07 [1] CRAN (R 4.1.0)
## gprofiler2 * 0.2.1 2021-08-23 [1] CRAN (R 4.1.0)
## gridExtra 2.3 2017-09-09 [1] CRAN (R 4.1.0)
## gtable 0.3.0 2019-03-25 [1] CRAN (R 4.1.0)
## haven 2.4.3 2021-08-04 [1] CRAN (R 4.1.0)
## highr 0.9 2021-04-16 [1] CRAN (R 4.1.0)
## hms 1.1.1 2021-09-26 [1] CRAN (R 4.1.0)
## htmltools 0.5.2 2021-08-25 [1] CRAN (R 4.1.0)
## htmlwidgets 1.5.4 2021-09-08 [1] CRAN (R 4.1.0)
## httpuv 1.6.3 2021-09-09 [1] CRAN (R 4.1.0)
## httr 1.4.2 2020-07-20 [1] CRAN (R 4.1.0)
## IRanges * 2.28.0 2021-10-26 [1] Bioconductor
## jquerylib 0.1.4 2021-04-26 [1] CRAN (R 4.1.0)
## jsonlite 1.7.2 2020-12-09 [1] CRAN (R 4.1.0)
## KEGGREST 1.34.0 2021-10-26 [1] Bioconductor
## knitr 1.36 2021-09-29 [1] CRAN (R 4.1.0)
## labeling 0.4.2 2020-10-20 [1] CRAN (R 4.1.0)
## later 1.3.0 2021-08-18 [1] CRAN (R 4.1.0)
## lattice 0.20-45 2021-09-22 [1] CRAN (R 4.1.0)
## lazyeval 0.2.2 2019-03-15 [1] CRAN (R 4.1.0)
## lifecycle 1.0.1 2021-09-24 [1] CRAN (R 4.1.0)
## limma 3.50.0 2021-10-26 [1] Bioconductor
## locfit 1.5-9.4 2020-03-25 [1] CRAN (R 4.1.0)
## lubridate 1.8.0 2021-10-07 [1] CRAN (R 4.1.0)
## magrittr 2.0.1 2020-11-17 [1] CRAN (R 4.1.0)
## MASS 7.3-54 2021-05-03 [1] CRAN (R 4.1.0)
## Matrix 1.3-4 2021-06-01 [1] CRAN (R 4.1.0)
## MatrixGenerics * 1.6.0 2021-10-26 [1] Bioconductor
## matrixStats * 0.61.0 2021-09-17 [1] CRAN (R 4.1.0)
## memoise 2.0.0 2021-01-26 [1] CRAN (R 4.1.0)
## mime 0.12 2021-09-28 [1] CRAN (R 4.1.0)
## modelr 0.1.8 2020-05-19 [1] CRAN (R 4.1.0)
## munsell 0.5.0 2018-06-12 [1] CRAN (R 4.1.0)
## mvtnorm 1.1-3 2021-10-08 [1] CRAN (R 4.1.0)
## numDeriv 2016.8-1.1 2019-06-06 [1] CRAN (R 4.1.0)
## pheatmap * 1.0.12 2019-01-04 [1] CRAN (R 4.1.0)
## pillar 1.6.4 2021-10-18 [1] CRAN (R 4.1.0)
## pkgbuild 1.2.0 2020-12-15 [1] CRAN (R 4.1.0)
## pkgconfig 2.0.3 2019-09-22 [1] CRAN (R 4.1.0)
## pkgload 1.2.3 2021-10-13 [1] CRAN (R 4.1.0)
## plotly 4.10.0 2021-10-09 [1] CRAN (R 4.1.0)
## plyr 1.8.6 2020-03-03 [1] CRAN (R 4.1.0)
## png 0.1-7 2013-12-03 [1] CRAN (R 4.1.0)
## prettyunits 1.1.1 2020-01-24 [1] CRAN (R 4.1.0)
## processx 3.5.2 2021-04-30 [1] CRAN (R 4.1.0)
## promises 1.2.0.1 2021-02-11 [1] CRAN (R 4.1.0)
## ps 1.6.0 2021-02-28 [1] CRAN (R 4.1.0)
## purrr * 0.3.4 2020-04-17 [1] CRAN (R 4.1.0)
## R6 2.5.1 2021-08-19 [1] CRAN (R 4.1.0)
## RColorBrewer * 1.1-2 2014-12-07 [1] CRAN (R 4.1.0)
## Rcpp 1.0.7 2021-07-07 [1] CRAN (R 4.1.0)
## RCurl 1.98-1.5 2021-09-17 [1] CRAN (R 4.1.0)
## readr * 2.1.0 2021-11-11 [1] CRAN (R 4.1.0)
## readxl 1.3.1 2019-03-13 [1] CRAN (R 4.1.0)
## remotes 2.4.1 2021-09-29 [1] CRAN (R 4.1.0)
## reprex 2.0.1 2021-08-05 [1] CRAN (R 4.1.0)
## rlang 0.4.12 2021-10-18 [1] CRAN (R 4.1.0)
## rmarkdown 2.11 2021-09-14 [1] CRAN (R 4.1.0)
## rprojroot 2.0.2 2020-11-15 [1] CRAN (R 4.1.0)
## RSQLite 2.2.8 2021-08-21 [1] CRAN (R 4.1.0)
## rstudioapi 0.13 2020-11-12 [1] CRAN (R 4.1.0)
## rvest 1.0.2 2021-10-16 [1] CRAN (R 4.1.0)
## S4Vectors * 0.32.2 2021-11-07 [1] Bioconductor
## sass 0.4.0 2021-05-12 [1] CRAN (R 4.1.0)
## scales 1.1.1 2020-05-11 [1] CRAN (R 4.1.0)
## sessioninfo 1.2.1 2021-11-02 [1] CRAN (R 4.1.0)
## shiny 1.7.1 2021-10-02 [1] CRAN (R 4.1.0)
## stringi 1.7.5 2021-10-04 [1] CRAN (R 4.1.0)
## stringr * 1.4.0 2019-02-10 [1] CRAN (R 4.1.0)
## SummarizedExperiment * 1.24.0 2021-10-26 [1] Bioconductor
## survival 3.2-13 2021-08-24 [1] CRAN (R 4.1.0)
## testthat 3.1.0 2021-10-04 [1] CRAN (R 4.1.0)
## tibble * 3.1.6 2021-11-07 [1] CRAN (R 4.1.0)
## tidyr * 1.1.4 2021-09-27 [1] CRAN (R 4.1.0)
## tidyselect 1.1.1 2021-04-30 [1] CRAN (R 4.1.0)
## tidyverse * 1.3.1 2021-04-15 [1] CRAN (R 4.1.0)
## tzdb 0.2.0 2021-10-27 [1] CRAN (R 4.1.0)
## usethis 2.1.3 2021-10-27 [1] CRAN (R 4.1.0)
## utf8 1.2.2 2021-07-24 [1] CRAN (R 4.1.0)
## vctrs 0.3.8 2021-04-29 [1] CRAN (R 4.1.0)
## viridisLite 0.4.0 2021-04-13 [1] CRAN (R 4.1.0)
## withr 2.4.2 2021-04-18 [1] CRAN (R 4.1.0)
## xfun 0.28 2021-11-04 [1] CRAN (R 4.1.0)
## XML 3.99-0.8 2021-09-17 [1] CRAN (R 4.1.0)
## xml2 1.3.2 2020-04-23 [1] CRAN (R 4.1.0)
## xtable 1.8-4 2019-04-21 [1] CRAN (R 4.1.0)
## XVector 0.34.0 2021-10-26 [1] Bioconductor
## yaml 2.2.1 2020-02-01 [1] CRAN (R 4.1.0)
## zlibbioc 1.40.0 2021-10-26 [1] Bioconductor
##
## [1] /Library/Frameworks/R.framework/Versions/4.1/Resources/library
##
## ──────────────────────────────────────────────────────────────────────────────